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Abstract 

We present first-principfes kinetic Monte Carlo (lp-kMC) simulations addressing the CO ox- 
idation reaction at Pd(100) for gas-phase conditions ranging from ultra-high vacuum (UHV) to 
ambient pressures and elevated temperatures. For the latter technologically relevant regime there 
is a long-standing debate regarding the nature of the active surface. The pristine metallic surface, 
an ultra-thin (\/5 x v / 5)-R27°PdO(101) surface oxide, and thicker oxide layers have each been sug- 
gested as the active state. We investigate these hypotheses with lp-kMC simulations focusing on 
either the Pd(100) surface or the PdO(lOl) surface oxide and intriguingly obtain a range of (T,p)- 
conditions where both terminations appear metastable. The predicted bistability regime nicely ties 
in with oscillatory behavior reported experimentally by Hendriksen and coworkers [Catal. Today 
105, 234 (2005)]. Within this regime we find that both surface terminations exhibit very similar 
intrinsic reactivity, which puts doubts on attempts to assign the catalytic function to just one 
active state. 
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I. INTRODUCTION 



Common use of Pd for catalytic exhaust gas purification 1 has motivated frequent studies 
on Pd single- crystal model catalysts to elucidate the underlying molecular-level mechanisms. 
Notwithstanding, despite significant efforts the surface structure and composition at near- 
ambient conditions remains unclear. This is largely due to inherent difficulties in achieving 
atomic-scale information/resolution in this technologically relevant regime of near-ambient 
pressures and at temperatures of 300-800 K. Generally, the outcome of such experiments 
seems strongly dependent on the preparation conditions and experimental setup. Specifi- 
cally, there exists a longstanding and controversial debate on whether reactivity is due to 
the pristine metal surface 2 , a surface oxide film 3 , or even a thicker bulk-like oxide overlayer. 4 

Under ultra-high vacuum (UHV) conditions, different oxidation stages of Pd(100) 5 have 
been thoroughly characterized both experimentally and theoretically. By low-energy electron 
diffraction (LEED) Chang and Thiel first identified five distinct ordered structures before 
the onset of bulk oxide formation: a p(2 x 2) oxygen adlayer, a c(2 x 2) adlayer, a (5 x 5), 
and a (y/E x y/E)R27° (for brevity henceforth denoted as \/o) reconstruction. 6 Todorova et 
al. 7 established that the latter \/5 structure corresponds to a single-layer of PdO(lOl) on 
top of Pd(100). 

At low pressures several experimental studies concluded that formation of this surface 
oxide is accompanied with a low catalytic activity. 8 ' 9 This was primarily attributed to the 
low CO binding energy at the surface oxide 5 , which Gao et al. estimated to be around 
0.5-0.6 eV. 10 At elevated pressures the situation is less clear. Lundgren et al. 11 compared 
structural information from in situ surface x-ray diffraction measurements in a pure oxygen 
environment to a systematic ab initio surface phase diagram 12 ' 13 , and concluded that the 
formation of bulk oxide on Pd(100) is kinetically severely hindered even at temperatures 
around 675 K. This suggests that at the (T, p)-conditions relevant for near-ambient CO 
oxidation the most likely surface terminations are either pristine Pd(100) with some coverage 
of O adsorbates, or a v5 monolayer. On top of this, reactor scanning tunneling microscopy 
(STM) measurements by Hendriksen et al. in fact suggested bistability in this near-ambient 
regime. 14-16 In the understanding of the work by Lundgren et al. this result refers to (T,p) 
conditions where both O@Pd(100) and the \fh surface oxide are metastable. 

Here, we investigate this hypothesis with first-principles kinetic Monte Carlo (lp-kMC) 
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simulations that focus on either of the two surface states, i.e. either O@Pd(100) or the 
\/5 surface oxide, thereby extending previous lp-kMC work that had focused exclusively on 
the \/5 phase 12 ' 17 . We indeed find a range of (T, p)-conditions where both surface models 
appear metastable. Moreover, in this bistability regime, both surface terminations yield 
roughly similar turnover frequencies (TOFs) for near-ambient pressures. In contrast and 
consistent with experiment we find the Pd(100) model to be the much more reactive one in 
the UHV regime. This puts experiments 2 ' 3 into perspective that claim one surface state to 
be the active one, irrespective of the specific gas-phase conditions. 



II. METHODS 

A. lp-kMC simulations 

Our central goal is to describe the range of gas-phase conditions in which the surface 
termination of the Pd(100) model catalyst is characterized by pristine Pd(100) or by the a/5 
surface oxide. In order to also account for stabilization effects through the ongoing reaction 
kinetics this requires analyzing the steady-state coverages of CO and O in microkinetic mod- 
els of these two surface terminations as a function of external feed conditions (T, PcO)Po 2 )- 
For this, two factors indicate that standard microkinetic modeling in terms of mean-field 
rate equations will not be sufficient: On pristine Pd(100) lateral interactions between ad- 
sorbed O and CO are known to be rather strong 18-21 , while on the y/E surface oxide it is 
its essentially one-dimensional geometric trench structure that will lead to adlayer inhomo- 
geneities. We thus opt for lp-kMC simulations as presently only technique that provides 
the desired microkinetic information while fully accounting for the correlations, fluctuations, 
and explicit spatial distributions of the chemicals at the catalyst surface. 22-24 

In lp-kMC the time evolution of the system is coarse-grained to the discrete rare-event 
dynamics. Relying on a Markov approximation rejection-free lp-kMC algorithms thus gen- 
erate state-to-state trajectories that in their average yield the probability-density function 
Pi{t) to find the system at time t in state i representing the corresponding potential en- 
ergy surface (PES) basin i. The propagation of this probability density is governed by the 
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Markovian master equation, 

dPi(t) 



(1) 



dt 



where the sums run over all system states j, and kij is the rate constant for going from 
state i to j. The central first-principles ingredients required for the lp-kMC simulation 
are thus the individual rate constants of the considered elementary processes (adsorption, 
desorption, diffusion, reaction). To keep this input tractable, i.e. to arrive at a finite 
number of inequivalent processes and corresponding rate constants, lp-kMC simulations are 
commonly performed on lattice models. In the following sections we will first provide the 
working equations to determine the first-principles rate constants and then detail the 
specific lattice models employed in the present work for CO oxidation at pristine Pd(100) 
and at the \/5 surface oxide. 

For given environmental (gas-phase) conditions (which specify the adsorption rate con- 
stants), the output of lp-kMC simulations are then the detailed surface composition and 
occurrence of each individual elementary process at any time. Since the latter comprises the 
surface reaction events, this also gives the catalytic activity in form of products per surface 
area and time (i.e. TOFs), either time-resolved, e.g. during induction, or time-averaged 
during steady-state operation. 

B. First-principles rate constants 

The first-principles rate constants in this work are evaluated following the approach put 
forward by Reuter and Scheffler 25 . In brief, this approach relies on (harmonic) transition 
state theory (TST) for bound to bound processes like diffusion, and kinetic gas theory 
together with detailed balance to calculate adsorption and desorption rate constants. As 
the approach and its derivation have been detailed before, we here restrict ourselves to the 
presentation of the working equations for self-containment. 

The adsorption rate constant for species % is given by the rate with which these particles 
impinge on the unit-cell surface area A uc and the local sticking coefficient S sti i(T), which 
gives the fraction of the impinging particles that actually stick to a given free site st at 
temperature T 




(2) 
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Here, ks is the Boltzmann constant, p% the partial pressure of species i, and the particle 
mass. In unactivated adsorption events the local sticking coefficient merely accounts for the 
number of inequivalent sites in the surface unit-cell (vide infra). For activated adsorption 
events or Eley-Rideal (ER) type CO oxidation events, it is additionally governed by the 
adsorption or reaction barrier AE i st j, respectively, 

where A st) i is a geometrical factor reflecting the relative share with which molecules impinge 
on the different inequivalent surface sites st. 

Desorption of a particle adsorbed on a surface site st is modeled as the time reversed 
process of adsorption, and its rate constant thus has to fulfill detailed balance or microscopic 
reversibility. 

k£AT,Pi) (AG st>i (T, Pi ) 



exp 
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where AG stt i(T,Pi) is the difference in Gibbs free energy between the particle adsorbed at 
the surface state and in the gas phase. This is approximated by the difference between 
the gas-phase chemical potential fi gas ,i{T, Pi) 26 and the binding energy of the particle in the 
adsorbed state E]^f. Following the procedure detailed in Ref. 26 we interpolate tabulated 
values 27 to determine the gas-phase chemical potentials at any gas-phase condition. 

The diffusion rate constant of an adsorbate from one surface site st to another site st' is 
approximated as 

«fi.(W¥W(-^? s i . (5) 



h J \ ksT 

where AE^ stli is the diffusion energy barrier. Langmuir-Hinshelwood (LH) type CO oxida- 
tion reactions are described with an equivalent expression containing the reaction barrier. 
They, as well as the ER reactions, are generally treated as associative desorption events 
though, i.e. the formed CO2 immediately desorbs from the surface. 

Within this approach the required first-principles input to determine the rate constants 
reduces to the binding energies of the species at the surface sites (for the desorption rate 
constant), as well as to their diffusion and reaction barriers. The latter require the determi- 
nation of the corresponding transition states, which is at present commonly achieved with 



standard transition state search algorithms like nudged-elastic band 28 ' 29 . Unfortunately, we 
found the low PES corrugation for lateral shifts of the \/5 overlayer 7,30 to severely affect the 
performance and reliability of such methods. We therefore instead identified the transition 
states through scans along suitable reaction coordinates. For diffusion barriers the specific 
reaction coordinate employed was the lateral coordinate along a straight line connecting the 
known initial and final state. For CO oxidation reactions the distance between the C atom 
of the CO molecule and the O adsorbate was employed. We carefully checked for hysteresis 
effects and estimate the uncertainty in the determined barriers to be of the order of ±0.1 eV. 
The effect of this energetic uncertainty on the presented lp-kMC results will be critically 
discussed below and is not found to affect any of our conclusions. 

All required total energies were obtained from density-functional theory (DFT) with the 
generalized gradient approximation functional by Perdew, Becke, and Ernzerhof (PBE) 31 
to treat electronic exchange and correlation. Using the plane-wave code CASTEP 32 with 
standard library ultrasoft pseudopotentials systematic convergence tests showed that the 
quantities of interest (binding energies, diffusion and reaction barriers) are converged to 
within 30 meV at the employed energy cut-off of 400 eV and fc-point density of 0.4 A -1 . 
The surfaces were modeled in periodic supercell geometries, containing vacuum separations 
of more than 10 A. For the a/5 surface oxide the calculations were done within a (1 x 1) 
surface unit-cell of surface oxide on top of four layers of Pd(100). For pristine Pd(100) 
we employed a (2 x 2) surface unit-cell and also four metal layer slabs. In both cases, all 
geometries were fully relaxed to residual forces below 50 meV/A, while keeping the bottom 
two slab layers fixed at their bulk positions. 

C. lp-KMC lattice models 

1. Pd(100) 

In setting up the lattice model for Pd(100) we exploit the limited coverage range, for 
which the lp-kMC simulations need to provide a faithful representation. Detailed experi- 
mental work 6 indicates the formation of surface oxides at around a critical O coverage of 
0.5 monolayer (ML, defined as the ratio of O to top layer Pd atoms). As the purpose of 
this work is to derive the gas-phase conditions under which no surface oxide is formed, the 
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FIG. 1. (Color online) Top view of Pd(100) and the y/b surface oxide, illustrating the employed 
lattice models. Hollow sites are depicted by yellow circles and bridge sites by green rectangles (see 
text). Pd(100) top- layer atoms are represented by large dark green spheres, Pd atoms in the y/5 
oxide by large light green spheres, and O atoms by small red spheres. 

developed model only targets the coverage range 9 < 0.5 ML and the elementary processes 
occurring in it. Similarly, CO coverages of more than 0.5 ML have only been characterized 
for strongly overstoichiometric CO pressures and low temperatures 19,33-35 . As this regime is 
outside the catalytic context, only CO coverages below 0.5 ML are of interest here, too. 

In this coverage range O atoms bind preferably at the fourfold hollow sites at Pd(lOO) 18,36 , 
while CO binds at twofold bridge sites 13,37 . In consequence we set up the lattice model shown 
in Fig. 1, in which O can exclusively occupy hollow sites and CO bridge sites. Lateral inter- 
actions between oxygen adatoms have been systematically calculated by Zhang, Blum, and 
Reuter and were found to be strongly repulsive at nearest and moderately repulsive at next- 
nearest hollow- hollow distance 18 . For the targeted coverage regime and intermediate temper- 
atures up to 600 K, this is taken into account in the lp-kMC model by blocking adsorption 
or diffusion events into configurations that would result in O atoms in nearest-neighbor po- 
sitions. For dissociative adsorption this implies a required motif of eight empty sites (two 
next-nearest neighbor sites for the actual adsorption and their immediately adjacent six 
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nearest-neigbor sites) known as the 8-site rule 21,38 40 that is required for an adsorption event 
to take place. 

Similarly repulsive interactions between adsorbed CO molecules are indicated by the 
experimentally characterized superstructures 19 ' 33-35 and DFT calculations 12 . Accordingly 
and accounting for the shorter bridge-bridge distances, we also block CO adsorption and 
diffusion events into configurations resulting in CO molecules up to next nearest-neighbor 
bridge-bridge distance. We note, however, that the next nearest-neighbor lateral repulsion 
is not too strong 40 , and might hence be overestimated by the simple site-blocking in our 
default model, henceforth denoted as Pd(100)-2NN. We assess this below by additionally 
considering a model, in which CO-CO coadsorption is blocked only at nearest bridge-bridge 
distance (Pd(lOO)-lNN) and find this not to affect our conclusions towards bistability. 

Finally, the absence of mixed O-CO adsorbate structures at low temperatures 33 suggests 
O-CO lateral interactions to also be repulsive. Unfortunately, no systematic first-principles 
studies characterizing these interactions have been performed to date. In this situation we 
note that the nearest-neighbor hollow-bridge distance at Pd(100) is smaller than the O-CO 
distance in the CO oxidation transition states described below. O-CO lateral interactions 
are thus modeled by blocking adsorption and diffusion events into configurations involving 
O and CO at such close distance. 

Within this lattice model and site-blocking rules we consider all non-concerted O and 
CO adsorption, desorption, as well as nearest-neighbor site diffusion and reaction events. 
Oxygen adsorption is modeled as dissociative process into two vacant next-nearest hollow 
sites. CO adsorption occurs unimolecularly into one free bridge site. Both adsorption 
processes are treated as non-activated, i.e. using a local sticking coefficient of 1 for the 
dissociative adsorption of 2 and 0.5 for CO, the latter accounting for the two bridge sites 
per surface unit-cell. For oxygen adsorption in the low-coverage range this is supported 
by explicit DFT sticking coefficient calculations 41 , while for CO we verified that lifting the 
CO molecule vertically up from the bridge adsorption site yields a path without adsorption 
barrier. This together with the observation that a CO lowered above an adsorbed oxygen 
atom is efficiently steered into a neighboring vacant bridge site, supports the assumption 
that all CO molecules impinging in the vicinity of a vacant bridge site will stick. 

Neglecting any lateral interactions beyond the short-range blocking rules, the required 
O and CO binding energies, diffusion and reaction barriers are independent of the local 
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TABLE I. Summary of DFT binding energies, diffusion and reaction barriers (LH and ER) used 
in the Pd(100) and y/5 lp-kMC models. All values are in eV. 



PdOy^ desorption barriers: 






E° on-site energy, V nearest-neighbor 


lateral interaction 11 




J7}Q 

O, bridge O , hollow 

-0.51 a -1.95 a 


CO, bridge 

-1.40 a 


^CO, hollow 

-1.92 a 


Vb-O.br-br ^0-0,hol-hol 

0.08 a 0.07 a 


Vb-O.br-hol 

0.08 a 




VcO-CO,br-br ^CO-CO ,hol-hol 

0.08 a 0.13 a 


VcO-CO,br-hol 

0.14 a 




Vb-CO,br-br Vb- CO, hoi -hoi 

0.06 a 0.11 a 


Vb-CO,br-hol 

0.13 a 


Vb-CO,hol-br 

0.12 a 


PdO y/E diffusion barriers 






CO,br-j.br CO.hol^hol 

0.4 a 0.6 a 


CO,br--rhol 

0.3 a 




O.br^br O.hol^hol 
1.2 a 1.4 a 


0,br-»hol 

0.1 a 




PdOy5 reaction barriers (LH) 






-^Obr.CObr ^Ohol.COhol 

1.0 a 1.6 a 


^Obr.COhol 

0.5 a 


^Ohol.CObr 

0.9 a 


PdO\/5 reaction barriers (ER) 






^Oholl Sohol2 

0.8 0.5 


^Obrl 

0.0 


^Obr2 

0.0 


Pd(100) desorption barriers 






E° 

O , hollow 

-1.25 


CO, bridge 

-1.93 






Pd(100) diffusion barriers 






CO,br-i>br 0,hol->hol 

0.14 0.28 






Pd(100) reaction barriers 






^CObr.Ohol 

0.9 



Rogal et at, PRB 77, 155410-12, (2008) 



adsorbate environment. The corresponding values computed within our DFT setup are 
summarized in Table I and are generally in good agreement with previous computations 7 ' 12 ' 18 . 
The obtained LH reaction barrier of 0.9 eV in particular is in good agreement with previous 
DFT calculations at varying coverages (0.76-1.05 eV) 42-44 . Attempts to calculate an ER type 
reaction path over a (2 x 2) oxygen adlayer at Pd(100) showed that CO and O coadsorption 
is energetically more favorable and thus an ER reaction is excluded from the Pd(100) model. 
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As apparent from Table I the O and CO diffusion barriers are very low. This high mobility 
severely limits the numerical efficiency of the lp-kMC simulations, which are completely 
dominated by frequent executions of diffusion events at minute time increments. In order to 
speed up the lp-kMC simulations we thus artifically raised the diffusion barriers by as much 
as 0.5 eV. As diffusion is then still fast enough to achieve full equilibration of the adlayer in 
the temperature range of interest to this study, we found this not to affect the computed 
coverages or TOFs at all. 



2. y/5 surface oxide 

For the \/5 surface oxide we employ the lp-kMC model established and detailed by Rogal, 
Reuter, and Schemer 12,13 . In brief, this model considers two non-equivalent sites named 
bridge and hollow as depicted on the right of Fig. 1. Since there is no process involving sites 
from adjacent bridge-hollow trenches, the lattice may be viewed as quasi one-dimensional. 
In analogy to the just described Pd(100) model, the elementary process list consists of all 
non-concerted adsorption, desorption, diffusion and LH reaction processes involving these 
sites. Dissociative O2 adsorption requires two neighboring sites, and is only hindered by 
a sizable adsorption barrier of 1.9 eV in the case of adsorption into two bridge sites. CO 
adsorption is unimolecular and not hindered by adsorption barriers. 

As only modification of the Rogal model we additionally consider an ER reaction mech- 
anism, as recently suggested by Hirvi et al. for bulk PdO(lOl) 4 . Reaction path calculations 
vertically impinging a CO molecule over O atoms adsorbed in the different surface sites in- 
deed also yield rather low reaction barriers over the v5 surface oxide, namely about 0.7eV 
over hollow and essentially zero over bridge. The latter rather astonishing result has very 
little consequences for the surface oxide stability and the catalytic activity in near stoichio- 
metric feeds though. Under corresponding gas-phase conditions the O coverage of bridge 
sites is negligible and ER reaction events involving O in bridge sites correspondinlgy made 
no significant contribution to the overall TOF in the simulations. 
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D. kMC simulation setup 



Both lp-kMC models were implemented using the kmos framework 45 . All simulations 
were performed in simulation cells containing (20 x 20) unit cells and using periodic boundary 
conditions. Systematic checks showed that the quantities of interest here, i.e. the average 
steady-state coverages and TOFs, are perfectly converged at these cell sizes. For defined 
gas-phase feed conditions (T, Pco,Po 2 ) the lp-kMC simulations eventually reach a steady 
state, with constant TOF and average surface coverages Qi jS t of species % on site st, 

Vi,st - v- a , > W 

where n denotes the value at the nth kMC step after steady state has been reached. Starting 
from different initial adsorbate distributions (clean surface, O-poisoned surface and CO- 
poisoned surface) and using different random number seeds we validated that this steady 
state is well defined, i.e. we never observed multiple steady-states in this system. In the 
kinetic phase diagrams shown in Fig. 2 below the obtained results are summarized in form 
of contour plots. All contour plots were interpolated using radial basis functions in order to 
reduce numerical noise. 



III. RESULTS 

A. Kinetic phase diagrams at 600 K 

Predicted oxygen and CO coverages on the v5, as well as on the Pd(100) surface are 
shown in Fig. 2 for a temperature of 600 K. For the Pd(100) surface the coverages refer to 
CO in bridge sites and O in hollow sites, as these are the only available sites for these species 
in the employed model, respectively. For the \/5 surface only the O and CO coverage in 
the upper hollow sites is shown, as this is the expected critical descriptor for the stability 
of the surface oxide layer 12 . For Pd(100) a coverage of 1 monolayer (ML) corresponds to 
one adsorbate (O or CO) per surface Pd atom. For the \/5 surface a coverage of 1 ML 
corresponds to a complete occupation of the upper hollow rows by the respective species. 

Both models yield the intuitively expected coverage variations with partial pressure. 
Starting with pure O2 gas conditions (going horizontally along the bottom of the panels in 
Fig. 2), both models yield the correct surface coverages in the respective stable regimes 
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FIG. 2. (Color online) lp-kMC computed O and CO coverage as a function of partial pressures for 
T = 600 K at Pd(100) (upper panels) and at the y/b surface oxide (lower panels). At Pd(100) the 
coverage refers to O in hollow sites and CO in bridge sites, at the y/5 the coverage refers to O and 
CO in upper hollow sites (see text). The thick black line in the left panels indicates the boundary 
where each phase is expected to be stable: For Pd(100) this corresponds to the p (2 x 2) oxygen 
adlayer, for PdO\/5 this is a coverage > 0.9 ML. 



by construction: For Pd(100) increasing O2 content in the gas phase leads to a gradual 
coverage increase starting from the clean surface. For the v5 the full ML O coverage in 
the bottom right part of the panel in turn reflects a fully intact layer of the surface oxide 
in which all upper hollow sites are fully covered by oxygen. Increasing CO pressure also 
influences the surface coverages as intuitively expected: An increasing CO population at the 
surface defers the stabilization of surface O at Pd(100) to higher O2 pressures, while it leads 
to a quicker depletion of upper hollow O atoms at the \/5 at decreasing O2 pressures. For 
highly overstoichiometric CO pressures (upper horizontal line in the panels) the Pd(100) 
model displays CO coverages of 0.5 ML. This is less than the highest CO concentration 
characterized experimentally, namely 0.75 ML 35 . This deviation results from the employed 
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FIG. 3. (Color online) Deduced region of bistability, i.e. gas-phase conditions where the lp-kMC 
models would predict the simultaneous stability of pristine Pd(100) and the v5 surface oxide (see 
text). All diagram show the same range of gas phase chemical potentials. 

simplified lateral interaction model, which likely features too repulsive CO-CO interactions. 
We will scrutinize this with the less repulsive Pd-INN model below, but also note that 
these CO-rich gas-phase conditions (and dense CO adlayers) are not the focus of our present 
interest. The same holds for the incorrect limit of a purely CO-covered a/5 surface in the 
upper left part of the respective panel (where the surface oxide would in reality be reduced 
away), as well as for the c(2 x 2) 0.5 ML O coverage in the lower right part of the Pd(100) 
panel (where instead a surface or bulk oxide would be formed). 

The very transition between CO- and O-covered regimes requires closer inspection. In 
the case of the Pd(100) surface a CO-poisoned surface would be catalytically inactive, but 
certainly stable. This surface is expected to be stable at least until the oxygen concentration 
does not rise above ~ 0.25 ML as this corresponds to a coverage regime representative for 
the experimentally characterized p(2 x 2) overlayer 11 ' 46 . In Fig. 2 we therefore denote the 
stability boundary of the Pd(100) surface at 0.25 ML coverage, i.e. we would expect a stable 
Pd(100) surface for any gas-phase conditions to the upper left of this line. Note that due 
to the steep coverage rise in the transition region, this stability boundary would be barely 
affected on the scale of Fig. 2 if we had e.g. chosen 0.4 ML coverage as the stability criterion. 
In the same spirit Rogal et al. have used a coverage exceeding 0.9 ML as stability criterion 
for the a/5 surface oxide before 13 . The corresponding stability boundary for the surface 
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oxide is also drawn in Fig. 2, and we would expect the surface oxide to be stable anywhere 
to the bottom right of this line. Again, on the scale of Fig. 2 it would make little difference, 
if a stability criterion of e.g. 0.95 ML or 0.99 ML coverage had been used. With the two 
stability regions thus quite narrowly defined, the central and intriguing feature of Fig. 2 is 
that there is a finite range of partial pressures where both models are predicted to be stable. 
This is highlighted again in Fig. 3, where the resulting bistability region is marked in green. 

B. Variation with temperature 

The results presented in the preceding section were obtained for T = 600 K, a temperature 
that falls in the middle of the temperature range from ~ 300—800 K that is generally most 
relevant for CO oxidation catalysis. Not least to make contact with the dedicated reactor 
STM experiments performed by Hendriksen et al. at 408-443 K 15 , an important next step 
is to assess the variation of our findings, in particular the existence of a bistability region, 
with temperature. Correspondingly and following the same protocol as before, the two lower 
panels of Fig. 3 summarize our findings for T = 500 K and T = 400 K. Both panels show the 
same range of gas-phase chemical potentials as the upper panel for 600 K, allowing to ascribe 
all differences between the three panels directly to kinetic effects, that is deviations from 
thermodynamic scaling. Intriguingly, these deviations are rather small, i.e. the differences 
for the three temperatures amount primarily to the expected scaling of the pressure axes 25 . 
This leaves the central result in form of a bistability over a finite range of gas-phase pressures 
essentially untouched within the temperature range 400-600 K. Particularly intriguing is 
that certainly for T = 600 K the catalytically most relevant near-ambient conditions at 
stoichiometric feed fall right within the region of bistability. For the lower temperatures the 
bistability region shifts increasingly towards O-rich conditions, such that e.g. the conditions 
of the Hendriksen experiments fall just at the border of this region. 

IV. DISCUSSION 

A. Origin and robustness of bistability region 

An immediate concern with the observed bistability region is that it is rather narrow 
in (po 2 ,p co )-sp&ce. Considering the range of uncertainties underlying the lp-kMC models, 
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FIG. 4. (Color online) "Constrained" ab initio thermodynamic phase diagram of the Pd(100)/\/5 
system, considering only ordered structures consistent with the lp-kMC lateral interaction models 
(see text). The thick dashed line denotes the stability boundary between adsorption phases on 
Pd(100) (upper left part) and adsorption phases on the v5 surface oxide (bottom right part). 
Additionally shown as hatched area is the stability region of Pd(100) phases with O coverage 
below 0.5 ML. This is the thermodynamic equivalent to the bistability region identified in the lp- 
kMC simulations, if the low O coverage defers formation of the (thermodynamically preferred) 
surface oxide. 

prominently the semi-local DFT energetics and lateral interaction model, this raises doubts 
as to the robustness of this finding. Fortunately, in the present case an analysis of the 
atomic-scale reason behind the bistability shows that its actual existence emerges rather 
independently from these uncertainties, which is why we lead these two discussions jointly 
in this section. 

A first important step towards an understanding of the atomic-scale origin of the bista- 
bility comes from the observed almost perfect thermodynamic scaling of the location and 
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extent of the bistability region in (po 2 5Pco)-space in the temperature range 400-600 K, cf. 
Fig. 3. This suggests that the actual reaction kinetics, and the uncertainties in the concomi- 
tant reaction barriers, is not central to its existence. This view is confirmed by the fact that 
an equivalent bistability region is already obtained within a "constrained" ab initio ther- 
modynamics approach 47,48 , i.e. an approach that neglects the reaction kinetics completely. 
Figure 4 shows the corresponding phase diagram for the same range of chemical potentials 
also underlying the panels in Fig. 3. For a comparison to the lp-kMC results, exactly the 
same DFT energetics and only ordered structures consistent with the lp-kMC lateral inter- 
action models are used. This explains small differences with respect to the corresponding 
phase diagram published before by Rogal, Reuter, and Scheffler 12 , which e.g. considered 
additional experimentally characterized CO ordered phases (not of relevance for the present 
discussion). Denoted by the thick dashed line is the stability boundary between adsorption 
phases on the pristine Pd(100) surface (upper left part) and adsorption phases on the y/E 
surface oxide (lower right part). For all catalytically relevant gas-phase conditions, this 
boundary runs in fact between CO-covered Pd(100) and different adsorption phases on the 
\/5, i.e. with increasing O-content in the gas phase constrained thermodynamics predicts an 
abrupt phase transition from a CO-poisoned surface directly to the surface oxide. Knowing 
that a critical local coverage of about ~ 0.5 ML O is necessary to induce the formation of 
the surface oxide, it is therefore interesting to also include in the phase diagram the stability 
region of O-containing Pd(100) phases with less than 0.5 ML O coverage, i.e. phases where 
(coming from the pristine metal side) oxide formation would not yet start. Intriguingly, the 
resulting region shown in Fig. 4 extends over similar gas-phase conditions as the bistability 
region deduced from the lp-kMC simulations. 

As such we ascribe the bistability region to gas-phase conditions, where the v5 surface 
oxide is thermodynamically more stable, but where a too low O coverage at Pd(100) would 
not readily induce the formation of the oxide. When crossing the bistability region from O- 
rich to CO-rich gas-phase conditions the system will thus prevail in the (thermodynamically 
preferred) oxidized state, while in the opposite direction it will prevail in the metal state as 
kinetic limitations to stabilize enough oxygen at the surface prevent the formation of the 
surface oxide. With this understanding of the atomic-scale origin, it is primarily the adsor- 
bate binding energies and the lateral interaction model that are crucial to the robustness 
of the bistability region - as they govern the stabilization of oxygen at the surface. More 
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specifically, it is the relative binding energy differences at the metal and the \/5 that will 
primarily affect the extension of the bistability region, while the absolute binding energet- 
ics will rather shift its location in (po 2 , pco) _s P ace - We would expect uncertainties in the 
employed DFT exchange-correlation functional to rather affect the absolute binding ener- 
getics (i.e. systematic over- or underbinding) and only to a lesser extent the binding energy 
differences of O and CO at the two surfaces. This view is confirmed by systematic tests, 
in which we increased or decreased all binding energies by 0.2 eV and in both cases found 
only small differences in the extension of the bistability region at more pronounced changes 
in its position in (p 02 ,pco)spa-ce. Correspondingly, we do not expect that the predicted 
gas-phase conditions for the bistability are accurate to better than some orders of magni- 
tude in pressure. Its actual existence, however, should be very robust with respect to the 
uncertainties of present-day DFT functionals. 

This leaves as final important aspect the employed lateral interaction model on the 
Pd(100) surface. Already the comparison to the experimentally determined maximum CO 
coverage (0.75 ML 35 versus 0.5 ML in Fig. 2) indicated the site-blocking rules in the default 
Pd(100)-2NN model to be too repulsive. We assess the consequences for the bistability re- 
gion by determining this region with the less repulsive Pd(100)-1NN model, where now only 
configurations that would result in CO adsorbed in nearest-neighbor positions are excluded. 
The resulting lp-kM C bistability region at 600 K is shown in Fig. 5 and compared to the 
corresponding bistability region deduced before within the default Pd(100)-2NN model. Ob- 
viously the eased stabilization of CO at the Pd(100) surface (maximum coverage now 1 ML) 
blocks the O adsorption even more effectively and consequently leads to a much increased 
range of gas-phase conditions where both surface states are predicted to be stable. If we 
take the comparison to the maximum CO coverage from experiment as measure, we would 
expect the true CO-CO repulsion to be somewhere between the one represented within the 
Pd(100)-1NN and Pd(100)-2NN site blocking models. Correspondingly, we would estimate 
the true extension of the bistability region to be somewhere between those predicted within 
the two models. While thus not fully quantified in its extension, the actual existence of a 
bistability region per se is robust against the uncertainty in the employed lateral interaction 
model. 
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FIG. 5. (Color online) Comparison of the bistability regions arising from lp-kMC simulations 
with different lateral interaction models. Pd(100)-2NN is the default model employed sofar, where 
configurations are excluded that would result in CO adsorbed in nearest-neighbor and next-nearest 
neighbor positions (green solid area). Pd(100)-1NN mimics less repulsive interactions and only 
excludes configurations resulting in CO adsorbed in nearest-neighbor positions (hatched area). 

B. Comparison to experiment 

The finding of a finite bistability region in our calculations agrees nicely with the bista- 
bility and oscillations reported experimentally by Hendriksen and coworkers 14 ' 15 . Unfortu- 
nately, we can not directly compare the temperatures and pressures where this bistability 
is found. On the theoretical side this is due to the aforediscussed uncertainty with respect 
to the location of this bistability region in (po 2 , j9co)- s P ace caused by the approximate first- 
principles parameters and lateral interactions entering the kinetic models. On the experi- 
mental side this is due to the suspected non-negligible mass-transfer effects in the employed 
reactor STM setup 49-51 . Nevertheless, the deduced pressure ranges are intriguingly close, as 
is the extension of the bistability region. At the experimental T ~ 400 K and po 2 — 1 atm, 
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the measured width of the bistability region in CO pressure is 0.02 atm 15 , while it is ~ 0.5 atm 
in the simulations. Also, the finding that at constant po 2 , the bistability region shifts with 
increasing temperature to increasing pco agrees with the experimental observations 14 . As 
such our interpretation is that the here obtained bistability between metallic Pd(100) and 
the a/5 surface oxide provides an atomic-scale model for the oscillations of the Hendriksen 
experiments. This is close to the interpretation arrived at in the experimental study, with the 
slight modification that they suspected the metallic phase to correspond to predominantly 
O-covered Pd(100), whereas our simulations show that this is predominantly CO-covered 
Pd(100). 

The atomic-scale insight into the surface composition provided by our lp-kMC simula- 
tions also allows to further qualify the reaction mechanism. An intriguing observation made 
by Hendriksen and coworkers was that in contrast to the metal surface, the reaction rate 
on the oxide did not follow traditional LH kinetics. They proposed this as a signature of 
low oxide stability and a concomitant Mars-van-Krevelen type mechanism, in which the 
oxide is continually consumed and reformed. Figure 6 shows that this can instead be ra- 
tionalized within a much simpler picture that does not involve the oxide stability. Shown 
are the surface coverages on the Pd(100) metal and \/5 surface oxide over the pco pressure 
range, where we obtain bistability at fixed T = 400 K and po 2 — 1 atm, i.e. for gas-phase 
conditions comparable to the Hendriksen experiments. At the Pd(100) surface, O and CO 
adsorb in inequivalent sites, hollow and bridge, respectively. Due to the strongly repulsive 
interactions, mimicked in our simulations by mutual site blocking, this nevertheless leads 
effectively to a competition for adsorption sites between the two adsorbates and in conse- 
quence to a LH-type kinetics. In contrast, on the a/5 surface oxide, O and CO also adsorb in 
inequivalent sites at these gas-phase conditions, namely O almost exclusively on the hollow 
sites and CO almost exclusively on the bridge sites, cf. Fig. 6. At the surface oxide these 
sites are located further away from each other and adsorbates even in nearest-neighbor sites 
experience only very small lateral interactions, cf. Table I. Consequently, the occupation of 
these two site types occurs almost independently of each other at varying gas-phase condi- 
tions and gives rise to a TOF proportional in CO pressure, cf. Fig. 6, exactly as observed 
in the experiment. Similarly, we find the TOF to vary only little when changing the oxygen 
partial pressure away from the conditions shown in Fig. 6 - again fully consistent with the 
experimental findings. 

19 



8 o.l- 




pco [atm] 



FIG. 6. (Color online) lp-kMC surface coverages and turn-over- frequencies (TOFs) at gas-phase 
conditions comparable to the reactor STM experiments by Hendriksen et al. 15 . Shown is the 
bistability region in pco at fixed T = 400 K and po 2 = 1 atm. Upper panel: O and CO coverage 
on the Pd(100) surface; middle panel: O and CO coverage on the \/5 surface oxide; lower panel: 
comparison of the intrinsic TOFs of the two system states. 
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The understanding that the actual reaction mechanism over the \/5 proceeds by adsorp- 
tion on the otherwise intact surface oxide suggests that the continued roughening of the 
surface observed in the experiments is a by-product and not essential to the catalytic ac- 
tivity. With respect to the total activity Fig. 6 shows that over the bistability region both 
surface states, Pd(100) and \/5, exhibit rather similar intrinsic TOFs to within a factor of 
three. While the mass-transfer limitations present in the reactor STM measurements pre- 
vent a direct comparison, this finding alone sets the controversial discussion concerning the 
active state of the surface into perspective. At low temperatures and UHV conditions, our 
lp-kMC models indeed yield a significantly higher activity of the Pd(100) surface. This arises 
predominantly from the comparatively weak CO binding at the \/E surface oxide and the 
concomitant limitations in stabilizing it at the surface. However, at technological gas-phase 
conditions, which in our simulations do fall within the bistability region, the TOF differences 
between the two models are not large enough and furthermore change sensitively with the 
detailed gas-phase conditions, cf. Fig. 6, to support a detailed discussion as to which state 
is the more active one. In this respect, the Pd(100) surface differs qualitatively from the 
equally prominently discussed Ru/Ru02 system 52,53 . On ruthenium the catalytic activity 
can be clearly attributed to the oxidized surface, while in the present system two competing 
surface states can equally contribute to the catalytic activity, which thus also explains why 
a clear signature of palladium (surface) oxide formation may in some circumstances be more 
difficult to find. 



V. CONCLUSIONS 

Detailed lp-kMC simulations on either the Pd(100) surface or the \/5 surface oxide yield 
a finite range of CO oxidation gas-phase conditions, where both surface states appear stable. 
This finding of a bistability region is robust against variations in the detailed criteria used 
to assess the stability of either state, as well as against the uncertainties arising from the 
approximate first-principles energetics and lateral interaction models. Notwithstanding the 
latter uncertainties do affect the position where in (T, p)-space the bistability region is found 
and we expect this to translate in uncertainties up to the order of 100 K and some orders of 
magnitude in pressure. 

The bistability region arises from limitations in stabilizing oxygen at the Pd(100) sur- 
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face, which extends the stability region of the metallic surface beyond that predicted by 
thermodynamics. In some respect this is therefore reminiscent of the generic Langmuir- 
Hinshelwood kinetics that lead to bistability between a CO-poisoned state that effectively 
blocks O2 adsorption and an O-rich state, i.e. the bistability of the reaction rate originates 
essentially from the inequivalence of the adsorption of carbon monoxide and oxygen. 54-58 . 
The difference is that in these classical models the two stable solutions were either supposed 
to be different adsorbate phases on the same substrate or did correspond to an active metal 
and an inactive oxide state. Instead, in the present system we have a CO-rich Pd(100) 
state that blocks oxygen adsorption and the surface oxide as O-rich state. Furthermore, the 
bistability does comprise technological conditions and conditions comparable to the dedi- 
cated reactor STM experiments performed by Hendriksen et a/. 15 ' 16 . Under these conditions, 
both surface terminations do show similar intrinsic activity, which sets preceding discussions 
with respect to the active state into perspective. 

Within the uncertainties of the exact location of the bistability region in (T, p)-space we 
suggest the here obtained bistability between metallic Pd(100) and the a/5 surface oxide 
as an atomic-scale model for the oscillations of the Hendriksen experiments. Likewise, our 
data lends further support to the view that oxide formation plays an important role in un- 
derstanding the catalytic activity of Pd catalysts. If it is even the oscillations themselves, 
catalytic activity would again be rationalized as a kinetic phase transition phenomenon, i.e. 
the catalyst surface being close to an instability 59-61 . At present our simulations performed 
separately on the two intact surface states, corresponding to ideal terraces of Pd(100) and 
the a/5 surface oxide, can not address this notion directly. This holds equally for the exper- 
imental interpretation that the continued roughening of the oxidic surface during reaction, 
as well as the formation of steps are crucial ingredients to the oscillatory behavior. Our 
ongoing work is directed to extend the lp-kMC capabilities towards the actual oxide for- 
mation and decomposition, which will then allow us to start scrutinizing these points with 
predictive-quality theory. 
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